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We  report  on  the  development  of  an  aerodynamic  force  model  for  a  flat-plate  wing 
undergoing  unsteady  motions  of  surging,  pitching,  and  the  combination  thereof  in  two- 
dimensional  flow.  This  investigation  aims  to  extend  the  conventional  quasi-steady  aero¬ 
dynamic  theory  to  account  for  the  influence  of  large-scale  vortices  generated  from  the 
massively  separated  flow  at  high  angles  of  attack.  A  data  set  compiled  from  direct  numer¬ 
ical  simulation  is  used  to  develop  the  force  model.  Particular  focus  is  placed  on  examining 
the  influence  of  large-amplitude  wing  motion  on  the  unsteady  aerodynamics  force.  This  is 
especially  important  since  the  resulting  circulatory  component  of  force  is  significant  in  size 
and  is  nonlinearly  related  to  the  flow  field  causing  deviation  from  conventional  aerodynamic 
theories.  The  present  force  model  is  constructed  to  predict  the  aerodynamic  force  on  the 
body  for  the  duration  of  the  wing  motion. 


I.  Introduction 

As  we  pursue  improvement  of  aircraft  maneuverability,  the  aircraft  is  expected  to  have  predictable  high 
aerodynamic  performance  at  high  angles  of  attack  and  to  operate  in  gusty  conditions.  Such  predictive 
capability  is  critical  for  unmanned  air  vehicles  that  are  designed  to  fly  in  urban  environment  where  the 
wake  generated  behind  buildings  and  atmospheric  fluctuations  pose  a  challenge  to  stable  flight  operations. 
In  order  to  achieve  desirable  flight  performance  at  these  conditions,  it  becomes  necessary  to  develop  an 
accurate  aerodynamic  force  model  for  unsteady  wing  motions  with  large  amplitudes.  Compared  to  past 
quasi-steady  type  force  theories,  it  is  important  to  account  for  the  massively  separated  region  that  forms 
behind  wings. 

The  objective  of  the  present  work  is  to  develop  an  unsteady  force  model  that  requires  only  wing  kinematics 
as  input  to  predict  the  normal  force  exerted  on  the  wing.  We  note  that  the  present  study  investigates  large 
amplitude  unsteady  wing  motions  that  departs  from  linear  aerodynamic  theories. As  we  consider  high 
angle  of  attack  maneuvers,  the  force  applied  on  the  body  is  nearly  normal  to  the  surface  of  the  wing  with 
L/D  «  1.  For  this  reason,  we  select  the  normal  force  to  be  modeled.  In  this  study,  we  limit  our  model  to  be 
what  we  refer  to  as  the  fast  model.  That  is,  we  do  not  consider  first  or  second-order  dynamics  in  the  model 
such  as  the  one  utilized  by  Goman  and  Khrabrov.’^ 

In  the  present  investigation,  we  base  the  development  of  the  force  model  on  the  understanding  gained 
from  past  studies  on  unsteady  aerodynamics  around  wings  undergoing  various  types  of  maneuvers.  The 
effect  of  unsteady  maneuvers  such  as  pitching,'^’ rotation,'^’ surging, 12  combinations  of 

these  motions^^  have  been  examined.  From  these  studies,  it  has  been  observed  that  the  rate  of  motion 
affects  the  aerodynamic  force  significantly.  It  also  has  been  observed  that  the  decreasing  aspect  ratio  of  the 
wing  has  an  effect  on  reducing  the  aerodynamic  force  with  more  influence  later  in  time  once  the  tip  vortices 
develop  significantly.®  On  the  other  hand,  the  change  in  unsteady  force  from  varied  Reynolds  number  can 
be  secondary  in  magnitude  for  Re  <  O(IO^). 
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For  this  study,  we  limit  the  motion  to  be  at  a  low  Reynolds  number  of  500  (based  on  terminal  velocity) 
composed  of  surging  and  pitching  that  is  defined  by  the  modified  Eldredge  function  which  does  not  introduce 
any  discontinuity  in  the  temporal  derivatives  of  the  motion  profile.  The  pitching  motion  is  taken  to  occur 
about  the  leading  edge.  The  flow  considered  here  is  two-dimensional  so  that  the  effect  of  aspect  ratio  can 
be  included  later  in  a  follow-up  study.  The  present  force  model  takes  the  wing  motion  as  the  sole  input 
to  predict  the  force  experienced  by  the  wing.  This  purpose  for  this  choice  of  input  is  to  allow  for  the 
construction  of  a  computationally  inexpensive  force  model  that  can  be  implemented  into  a  vehicle  control 
algorithm  in  a  simple  manner.  In  what  follows,  we  describe  the  flow  physics  around  moving  wings  and 
discuss  the  development  of  the  force  model,  which  accounts  for  the  effects  of  surging,  pitching,  and  incidence 
angle.  A  variety  of  motions  are  examined  and  comments  are  offered  on  the  validity  of  the  model. 

II.  Approach 

We  focus  on  unsteady  wing  motion  in  two-dimensional  flow  for  a  Reynolds  number.  Re  =  t/maxc/i^  =  500, 
where  Umax  is  the  maximum  final  velocity  attained  by  the  wing  at  the  end  of  the  motion  considered  in  this 
study.  This  Reynolds  number  is  chosen  to  keep  the  flow  laminar  and  to  allow  for  this  study  to  be  performed 
in  the  absence  of  turbulence.  It  should  be  noted  that  the  qualitative  behavior  of  the  flow  field  at  Reynolds 
number  below  10®  is  mostly  unchanged  for  unsteady  wing  motion,  which  allows  for  this  study  to  be  applicable 
to  higher  Reynolds  numbers.®  We  employ  the  immersed  boundary  projection  method^'^AS  ^  simulate  the 
flow  around  the  wing  and  the  associated  aerodynamic  force.  The  immersed  boundary  projection  method 
solves  for  the  flow  field  on  a  Cartesian  mesh  and  generates  a  body  of  arbitrary  geometry  with  a  set  of 
Lagrangian  points  where  boundary  forces  are  applied  to  enforce  the  no-slip  boundary  condition.  A  typical 
size  of  the  domain  used  in  this  study  is  {x,y)  G  [—10,22]  x  [—16,16],  with  the  leading  edge  of  the  wing 
position  positioned  at  the  origin.  The  spatial  coordinates  are  non-dimensionalized  by  the  chord  length,  c, 
and  the  temporal  units  are  reported  in  terms  of  the  convective  time.  The  forces  on  the  body  (Fx^Fy)  are 
non-dimensionalized  using  the  dynamic  pressure  based  on  the  final  maximum  translational  velocity  Umax 
and  the  chord  c,  defined  by  Cl  =  Fy/  {^pU^^^c)  and  Co  =  F^/  (^pU^ax*^))  where  p  is  the  density  of  the 
fluid. 

III.  Flow  Physics 


Unsteady  Maneuvers 

We  consider  the  motion  of  surging  (accelerating),  pitching,  and  the  combination  thereof  starting  from  rest. 
For  the  pitching  motion,  we  limit  our  discussion  to  pitching  about  the  leading  edge.  The  effect  of  pitching 
about  other  pivot-points  along  chord  of  the  wing  are  discussed  in  Ref.  16.  The  velocity  and  angle-of- 
attack  profiles  used  for  the  surging  and  pitching  motions,  respectively,  are  defined  by  the  modified  Eldredge 
function  for  a  smoothed  linear  ramp.  A  linear  ramp  in  velocity  corresponds  to  constant  acceleration  and  a 
quadratic  profile  for  the  position,  whereas  a  linear  ramp  in  angular  position  corresponds  to  constant  angular 
velocity.  The  constant  acceleration  portion  of  the  motion  used  in  this  study  is  similar  to  the  motions  studied 
by  Chen  et  al.,^'^  whose  velocity  profile  was  defined  by  a  power  law  of  the  form  U{t)  =  kt^ ,  with  /3  =  1. 

For  the  surging  motions,  the  velocity  profile  is  defined  with  respect  to  convective  time  {t  =  Ut' /c),  given 

by 


U{t) 


Umax  'cOSh(a(t  -  ti)) 

2aCs  _cosh(a(<  —  12)) 


max 

~Y~ 


(1) 


where  Umax  is  the  maximum  velocity  the  plate  achieves  at  the  end  of  the  motion,  a  is  a  smoothing  parameter 
for  the  linear-ramp  profile,  and  Cg  =  t2  —  U  is  the  duration  of  the  motion.  Once  the  velocity  of  the  plate  is 
known,  the  dimensional  (physical)  time  =  ct/U)  is  determined  from 


t'  =  c 


dr 

W)' 


(2) 


By  integrating  the  velocity  profile,  the  position  of  the  plate  is  determined  by  numerically  integrating  the 
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velocity  profile: 


s  = 


U  {T)dT. 


(3) 


Throughout  this  paper,  we  consider  the  position  of  the  plate  to  be  the  convective  variable  of  the  motion. 
This  allows  for  motions  to  be  described  as  occurring  over  a  certain  number  of  chord  lengths  in  which  the 
flat  plate  has  traveled.  For  the  force  modeling  in  the  following  sections,  the  models  in  consideration  take  the 
wing  motion  as  the  sole  input,  including  the  velocity  of  the  plate  and  its  derivatives  (s,  s, ...).  The  velocity 
of  the  plate  as  well  as  the  higher  derivatives  are  found  by  taking  consecutive  derivatives  of  the  position  of 
the  plate  (s)  with  respect  to  the  physical  time,  t' ,  i.e.,  s  =  ds/dt',  s  =  ...  .  The  position  (s),  velocity 

(s),  and  acceleration  (s)  of  the  plate  are  shown  plotted  against  convective  time  (t,  left)  and  physical  time 
(t',  right)  in  Figure  1. 


H=s)  t' 


Figure  1:  The  position  (black),  velocity  (blue),  and  acceleration  (red)  profiles  for  the  surge  motion  plotted  against  convective 
time  (t,  left)  and  physical  time  (£',  right).  For  the  present  motion,  Cs  =  t2  —ti  =  1,  which  corresponds  to  a  constant  acceleration 
over  one  chord  (c)  length  of  travel. 


For  the  pitching  motions,  the  angle  of  attack  of  the  flat  plate  with  respect  to  convective  time  is  given  by 


a{t) 


a^log  'cosh(Q(f  -  ti)) 

2acp  [cosh(a(f  — 12)) 


Q^max 

2 


(4) 


where  amax  is  the  maximum  angle  the  plate  reaches  with  respect  to  the  freestream  and  Cp  =  t2  —  ti  is 
the  duration  over  which  pitching  occurs.  The  derivatives  of  the  angular  position  which  are  needed  for  the 
modeling  are  also  taken  with  respect  to  physical  time,  t',  i.e.,  a  =  da/dt' ,  a  =  d?a/dt'^ , ... 


Vortex  Dynamics 

In  order  to  understand  how  the  aerodynamic  forces  are  generated  with  unsteady  plate  motions,  we  first  focus 
on  the  evolution  of  the  vortical  structures  around  the  flat  plate.  Figure  2  shows  the  evolution  of  the  flow 
field  around  the  flat-plate  wing  undergoing  a  pure  pitching  (left  column),  combined  pitch-surging  (center 
column),  and  pure  surging  (right  column)  motion.  The  rate  at  which  each  of  these  maneuvers  takes  place 
is  given  by  a  reduced  frequency  oi  K  =  acj^Uao  =  amaxc/2C/ooCp  =  tt/S,  which  corresponds  to  the  motion 
taking  place  over  1  chord  length  of  travel.  For  the  pure  pitching  motion,  the  plate  resides  at  0°  angle  of 
attack  for  10  convective  units  to  allow  for  the  initial  transients  of  the  simulation  to  die  out  and  for  the 
forces  exerted  on  the  plate  to  reach  a  steady  value.  For  the  pure  surging  motions,  the  plate  resides  in  a 
quiescent  flow  at  an  angle  of  attack  of  a  =  0(^0)  with  respect  to  the  acceleration  direction.  For  the  combined 
pitch-surge  motion,  the  flat  plate  begins  the  maneuver  from  rest  at  an  angle  of  attack  of  a  =  0°  with  respect 
to  the  direction  that  the  acceleration  occurs.  The  evolution  of  the  flow  field  is  visualized  using  contours  of 
constant  spanwise  vorticity  given  by  —25  <  Wz  <  25. 

Beginning  with  the  pure  pitching  motion  in  the  left  column,  the  plate  has  a  well  established  boundary 
layer  at  the  start  of  the  motion  which  rolls  up  and  detaches  at  the  trailing  edge  as  the  plate  increases  is  angle 
of  attack  which  then  begins  to  advects  downstream.  At  the  leading  edge  of  the  plate,  the  vortex  sheet  rolls 
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Figure  2:  Evolution  of  the  flow  field  around  a  flat-plate  wing  undergoing  a  pure  pitching  (left  column),  combined  pitch-surging 
(center  column),  and  pure  surging  (right  column)  motion  at  a  reduced  frequency  of  A"  =  tt/S  and  Re  =  500.  The  flow  around 
these  plates  is  visualized  using  contours  of  —25  <  ^2  ^  25. 
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up  into  a  compact  leading-edge  vortex  (LEV)  that  presents  a  low  pressure  region  on  the  top  surface  of  the 
wing.  This  low  pressure  region  enhances  the  lift  of  the  plate  as  it  continues  to  increase  its  angle  of  attack. 
By  the  end  of  the  motion  {t  >  1),  the  plate  is  at  an  angle  of  a„iax  =  45°  and  the  LEV  continues  to  grow 
over  the  surface  of  the  wing. 

Unlike  the  pure  pitching  motion,  there  is  no  well  established  boundary  layer  at  the  start  of  the  motion 
for  the  surge  and  combined  pitch-surge  maneuver  since  the  motion  of  the  plate  begins  from  rest.  This  results 
in  the  absence  of  trailing  edge  vortices  for  both  the  pure  sure  motion  and  the  combined  pitch-surge  motion. 
For  the  pure  surge  motion,  the  flow  detaches  from  the  leading  edge  at  the  beginning  of  the  motion  due  to 
the  plate’s  high  angle  of  attack  of  Omax  =  45°.  This  results  in  a  relatively  weak  LEV  that  forms  during  the 
motion,  which  is  seen  to  be  smeared  along  the  chord  of  the  wing  when  compared  to  the  LEV  that  develops 
for  both  the  pure  pitching  and  combined  pitch-surge  motions.  For  the  pitch-surge  motion,  the  plate  begins 
to  accelerate  from  a  =  0°,  which  results  in  the  flow  being  attached  during  the  early  potion  of  the  motion. 
As  the  plate  continues  to  accelerate  and  increase  its  angle  of  attack,  the  flow  field  resembles  that  of  the  pure 
pitching  motion.  By  the  end  of  the  motion,  the  LEV  appears  to  be  more  compact  and  located  closer  to  the 
wing  surface  with  a  greater  strength  compared  to  the  LEV  that  develops  for  the  pure  pitching  case,  which 
is  due  to  the  increased  circulation  due  to  the  combined  motion,  and  enhances  the  lift  force. 


Aerodynamic  Forces 

The  temporal  evolution  of  the  lift  and  drag  coefficients  for  the  three  motions  under  consideration  is  given  in 
Figure  3  for  both  the  present  simulations  {Re  =  500,  solid  curve)  and  companion  experimental  measurements 
performed  in  a  horizontal  water  tunnel  at  the  Air  Force  Research  Laboratory  at  Wright-Patterson  Air  Force 
Base^®(i?e  =  20,000,  dashed  curve).  As  the  plate  begins  to  increase  its  angle  of  attack  for  the  pure  pitching 
motion,  a  large  spike  in  the  lift  coefficient  is  observed,  which  results  from  the  initial  angular  acceleration  of 
the  plate.  As  the  plate  continues  to  increase  its  angle  of  attack,  the  lift  force  begins  to  rise  again  due  to  the 
development  of  the  LEV  over  the  top  surface  of  the  wing.  At  the  end  of  the  motion,  a  negative  spike  in  the 
lift  coefficient  is  observed  due  to  the  angular  deceleration  of  the  plate. 


t 


t 


Figure  3:  Comparison  between  the  experimental  {Re  =  20, 000,  dashed)  and  present  simulation  {Re  =  500,  solid)  force 
coefficients  for  a  flat  plate  wing  undergoing  a  pitching  (red  curve),  surging  (black  curve),  and  a  combined  pitch-surge  motion 
(blue  curve)  at  a  reduced  frequency  oi  K  =  tt/S,  which  corresponds  to  occurring  over  one  chord  length  of  travel. 


For  the  pure  surging  and  combined  pitch-surge  maneuvers,  the  forces  are  initially  zero  since  the  plates 
start  from  rest.  As  the  motion  begins,  the  lift  and  drag  coefficients  begin  to  increase.  This  increase  is  much 
more  rapid  for  the  combined  motion  due  to  the  very  tight,  compact  LEV  that  forms  above  the  wing.  Towards 
the  end  of  the  maneuver,  the  lift  coefficient  for  the  combined  motion  is  nearly  three  times  that  of  the  pure 
surging  motion.  Although  the  initial  evolutions  of  the  lift  and  drag  coefficients  for  these  three  motions  are 
vastly  different,  both  lift  and  drag  forces  on  the  plate  collapse  for  the  three  maneuvers  roughly  2  convective 
units  after  the  completion  of  the  motion. 

Similar  trends  between  the  lift  and  drag  coefficients  for  both  Reynolds  numbers  considered  are  observed. 
Although  there  is  two  orders  of  magnitude  difference  in  the  Reynolds  number  between  the  present  simulations 
and  the  experiments,  we  observe  agreement  between  the  forces  generated  for  these  high-amplitude  maneuvers. 
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This  observation  gives  us  confidence  that  the  force  modeling  performed  in  the  following  sections  for  the 
present  low-Reynolds  number  maneuvers  can  be  extended  to  higher  Reynolds  number  flows. 

IV.  Force  Model 


Pure  Surge  Motion 

First,  let  us  consider  the  wing  at  a  constant  angle  of  attack  of  a  =  45°  accelerating  from  rest.  The  wing  will 
accelerate  linearly  over  Cs  convective  units  before  reaching  a  constant  velocity  17max-  Based  on  prior  research 
on  unsteady  flows  around  moving  wings,  we  have  observed  that  the  normal  force  is  the  dominant  force  acting 
on  the  wing.  In  fact  for  most  conditions  considered  in  this  study,  the  lift  to  drag  ratio  is  approximately 
one.  For  this  reason,  we  will  model  the  normal  force  on  the  wing  in  the  present  study  for  all  motion  profiles 
considered. 

In  order  to  capture  the  normal  force,  we  have  considered  various  functional  forms  based  on  the  kinematics 
and  observed  that  the  following  expression  approximates  the  normal  force  well; 

Fat  «  cis(v^+ s),  (5) 

where  Ci  is  a  constant  to  be  selected  and  the  dot  represents  the  time  derivative.  The  right  hand  side  of  this 
equation  is  found  to  perform  fairly  well  for  a  wide  range  of  acceleration.  Based  on  this  observation,  we  have 
generalized  the  force  model  to  be 

Fjy  «  CiS^^  +  C2S^^S.  (6) 

We  test  this  model  for  1  <  <  8  and  tune  the  parameters  ci,  C2,  /3i,  and  /32  for  each  value  of  Cg  using  a 

least-squares  fit.  The  comparison  of  the  DNS  data  and  the  force  model  that  is  tuned  for  each  case  is  shown 
in  Figure  4.  It  can  be  seen  that  when  the  optimal  model  coefficients  are  used,  the  model  predicts  the  force 
so  well  that  they  are  almost  indistinguishable  during  acceleration.  The  error  during  acceleration  is  limited 
to  less  than  5%  of  the  maximum  normal  force  during  the  time  interval  of  interest.  While  comparing  the 
forces,  we  have  found  that  the  comparison  of  the  maximum  force  value  is  an  effective  indicator  to  evaluate 
the  accuracy  of  the  force  models.  In  what  follows,  we  present  the  comparison  of  the  force  model  with  the 
DNS  results  along  with  the  maximum  force  values  and  the  errors. 


t 


Figure  4:  The  force  model  based  on  Eq.  (6)  with  its  coefficients  tuned  for  each  surging  case  is  compared  with  the  DNS  results. 
The  left  figure  compares  the  force  model  (dashed)  against  the  DNS  results  (solid).  The  center  plot  compares  the  maximum 
normal  force  for  each  of  the  surging  case  (circle  =  DNS  and  square  =  model).  The  right  plot  is  the  error  between  the  model 
and  the  DNS  results  normalized  by  the  maximum  force  over  the  time  interval  considered. 


While  the  model  in  Figure  4  predicts  the  DNS  forces  quite  well,  it  would  be  more  desirable  to  find  the 
force  model  that  does  not  require  fitting  of  the  model  coefficients  for  each  and  every  case.  In  order  to  find 
the  trend  associated  with  the  effect  of  acceleration  on  the  normal  force  and  its  model,  let  us  curve  fit  the 
model  for  a  selected  surging  case  and  use  the  selected  model  parameters  (ci,  C2,  /3i,  and  /32)  for  all  other 
surging  cases.  We  show  the  comparison  made  for  models  tuned  at  Cg  =  2,4,  and  8  in  Figure  5.  It  can  be 
noticed  that  the  model  tuned  for  faster  surge  rates  do  not  predict  the  slower  surge  rate  cases  well.  The  large 
difference  can  be  seen  at  higher  Cg  for  the  model  tuned  at  Cg  =  2.  However,  the  model  tuned  for  slower  surge 
rates  predicts  the  faster  surge  rate  cases  well.  These  observations  suggest  that  some  of  the  model  parameters 
need  to  be  carefully  tuned  to  predict  the  normal  force  accurately. 
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optimized  for  Cg  =  8  optimized  for  Cg  =  4  optimized  for  c 


Figure  5:  The  force  model  based  on  Eq.  (6)  with  its  coefficients  tuned  only  for  the  indicated  surging  case  is  compared  with 
the  DNS  results.  See  Figure  4  for  legends  and  symbols. 
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Cl 

C2 


Figure  6:  The  variation  in  the  exponents  /3i  and  /32  and  the  coefficients  ci  and  C2  (left).  The  variation  in  coefficients  ci  and 
C2  with  fixed  exponents  (right).  Shown  are  the  linear  approximations  for  ci  and  C2. 


The  behavior  of  the  model  parameters  tuned  for  each  surge  rate  case  is  presented  in  Figure  6.  Except  for 
the  results  for  Cg  =  6,  the  values  of  the  exponents  Pi  and  P2  do  not  vary  significantly.  Hence,  we  approximate 
the  exponents  with  their  average  values 

Pi  — >■  avg/3i  =  1.89  and  P2  — >■  avg/32  =  0.138.  (7) 

For  the  coefficients  in  front  of  the  force  model  terms,  it  is  possible  to  use  the  mean  values  following  the 
approach  taken  for  the  exponents: 

Cl  — >■  avgci  =  2.03  and  C2  — >■  avgC2  =  2.23.  (8) 

We  can  also  consider  the  coefficients  ci  and  C2  to  be  linearly  related  to  Cg,  as  shown  in  Figure  6.  The  linear 
approximations  for  ci  and  C2  are 

Cl  — y  — 0.212cg  2.98  and  C2  — y  0.2636cg  F  1.04.  (9) 

These  coefficients  are  shown  with  the  exponents  pi  and  P2  fixed  to  the  mean  values. 

Next,  let  us  compare  the  force  model  performance  based  on  the  choices  proposed  above.  In  Figure  7,  we 
compare  the  models  based  on 

(i)  coefficients  tuned  for  the  surge  case  with  Cg  =  4, 

(ii)  using  mean  values  for  all  exponents  and  coefficients  (using  Eqs.  (7)  and  (8)),  and 

(iii)  using  the  mean  exponents  and  linear  model  for  the  coefficients  (based  on  Eqs.  (7)  and  (9)) 

Of  the  three  models  considered  here  and  shown  in  Figure  7,  we  find  that  setting  the  exponents  to  be  the 
mean  values  and  linearly  approximating  the  coefficients  accurately  predicts  the  unsteady  force  on  a  wing 
undergoing  surging  motion.  Thus,  we  use  the  last  force  model  in  what  follows.  To  summarize,  the  unsteady 
normal  force  model  for  the  surging  case  is 


(pure  surge)  Fm  «ciS^^  +  C2S^^s, 

where  ai  =  1.89,  a2  =  0.138, 

Cl  =  -0.212cg  +  2.98,  C2  =  0.2636cg  +  1.04. 


(10) 


Variation  in  the  Angle  of  Attack 

For  the  above  model,  we  have  only  considered  the  case  of  the  wing  being  at  an  angle  of  attack  of  45°.  In 
order  to  have  the  normal  force  prediction  be  applicable  to  any  angle  of  attack,  we  pre-multiply  the  above 
normal  force  model  by  a  factor  of  sin(a)/-\/2  so  that  the  model  intersects  the  results  at  45°.  This  also  ensures 
that  the  normal  force  becomes  zero  at  a  =  0°  and  achieves  the  maximum  value  at  90°,  respectively.  The 
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mean  exponents  only  all  mean  values  optimized  for  c 


a 


oi 

a 


a 


Figure  7:  Summary  of  the  model  comparison.  Shown  are  the  models  (top)  tuned  for  Cs  =  4,  (middle)  using  mean  values  for 
all  coefficients,  and  (bottom)  only  fixing  the  exponents  with  mean  values.  See  Figure  4  for  legends  and  symbols. 
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overall  results  for  the  example  of  Cs  =  4  for  various  angles  of  attack  are  presented  in  Figure  8.  It  can  be  seen 
that  this  simple  scaling  captures  the  normal  force  well  for  a  wide  range  of  incidence  angles: 


(11) 


Figure  8:  The  surging  case  results  shown  for  a  =  15°,  30°,  45°,  and  60°. 


Extension  to  Pitch-Surge  Motions 

Next,  we  extend  the  above  force  model  for  surging  to  a  combined  pitch-surge  motion  initiating  from  rest. 
We  first  consider  the  surge  and  pitching  motions  to  take  place  at  the  same  reduced  frequency  and  generalize 
it  for  motions  with  different  rates.  The  proposed  force  model  expressed  by  Eq.  (11)  includes  the  influence 
of  the  angle  of  attack  a.  Since  the  first  and  second  time  derivatives  of  the  incidence  angle  is  not  included  in 
this  model,  we  do  not  expect  to  accurately  predict  the  unsteady  force  especially  for  faster  maneuvers  where 
quasi-steady  assumptions  do  not  hold.  This  can  be  seen  from  Figure  9,  where  the  normal  force  from  fast 
motion  exhibit  large  difference  from  the  force  predicted  by  Eq.  (11). 
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Figure  9:  Naive  application  of  Eq.  (11)  to  combined  pitch-surge  motion  cases  (solid  lines:  DNS,  dashed  lines:  model). 

In  order  to  account  for  the  difference  between  the  aforementioned  force  model,  Eq.  (11),  and  the  forces 
generated  by  the  pitching  motion,  we  utilize  the  insights  obtained  from  research  on  pure  pitching  cases. 
That  is,  we  add  the  unsteady  forces  that  are  proportional  to  a  and  a.  Hence,  for  the  pitch-surge  force  model. 
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we  arrive  at 


(12) 


Fn  ~  sin(a)  [cis^^  +  C2S^^s]  +  dids^  +  dibi'^ , 

v2 


where  two  coefhcients  di  and  d2  need  to  be  empirically  determined.  In  Figure  10,  we  present  the  above 
model  with  di  and  d^  based  on  the  DNS  results.  The  other  coefficients  and  exponents  are  held  fixed  from 
the  pure  surging  cases.  This  proposed  pitch-surge  model  tuned  for  each  case  considered  in  Figure  10  agrees 
well  with  the  normal  force  obtained  from  simulation.  We  can  also  observe  that  coefficients  di  and  d2  can  be 
approximated  by  a  linear  function  and  the  mean  value,  respectively.  With  such  simplification,  we  arrive  at 
the  final  version  of  the  combined  pitch-surge  force  model  of 


(pitch-surge)  Fn  « 

where  /3i 
Cl 
di 


V2 


sin(a)  [cis^^  -I-  C2S^^s]  -I-  didP  +  d2dLS^ 


=  1.89,  152  =  0.138, 

=  -0.212cs  -h  2.98,  C2  =  0.2636cs  -f  1.04, 
=  0.834cp  -f  7.61,  d2  =  0.582. 


(13) 


Figure  10:  The  normal  force  prediction  on  the  left  (solid:  DNS  and  dashed:  model)  made  by  the  combined  pitch-surge  model 
of  Eq.  (12)  where  the  coefficients  di  and  d2  are  tuned  for  each  case,  as  shown  on  the  right. 
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Figure  11:  The  normal  force  prediction  based  on  the  combined  pitch-surge  model,  Eq.  (13),  with  the  coefficients  given  by 
mean  values  and  linear  function  of  Ca  (solid:  DNS  and  dashed:  model). 

This  final  form  of  the  model  captures  the  normal  force  on  a  wing  accurately  for  a  wide  range  of  pitch/surge 
rates,  as  shown  in  Figure  11.  Note  that  the  rates  shown  in  this  figure  (c^  =  Cp  =  1,  2, 4,  6, 8)  span  across  the 
rate  corresponding  to  the  universal  formal  time  of  vortices.  Hence  it  can  be  seen  that  this  proposed  model 
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can  approximate  the  force  very  well,  not  only  during  the  wake  vortex  formation  but  also  after  the  creation 
of  the  maximally  grown  wake  vortex. 

In  order  to  assess  how  the  model  performs  for  more  general  motions,  we  will  now  compare  the  present 
model  to  simulations  for  combined  pitch-surge  motions  where  the  pitch  and  surge  rates  are  different.  Figure 
12  shows  the  results  for  the  combined  pitch-surge  motions  over  a  range  of  pitch(cp  =  1,4,  and  6)  and  surge 

rates  (cg  =  1,4,  and  6).  The  vertical  lines  in  the  plots  highlight  regions  where  the  pitching  motion  ( - ) 

and  surging  motion  ( - )  are  completed.  The  opaque  portion  of  the  plots  indicate  regions  after  the  wing 

motion  has  completed  for  which  the  force  model  is  not  designed  to  predict.  The  deviation  in  force  after 
the  wing  maneuver  results  from  the  fact  that  the  only  input  to  the  model  is  the  wing  kinematics.  In  some 
cases,  however,  it  is  observed  that  after  the  motion  has  finished  the  model  continues  to  predict  the  force 
on  the  plate  reasonably  well.  This  is  primarily  due  to  the  inherent  dynamics  of  the  LEV  formation  that  is 
embedded  within  the  coefficients,  which  were  empirically  determined  to  match  the  DNS  results. 
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Figure  12:  The  normal  force  prediction  based  on  the  combined  pitch-surge  model,  Eq.  (13),  for  maneuvers  in  which  the 

pitch-rate  and  surge-rate  are  different.  The  vertical  lines  correspond  to  the  end  of  the  pitching  motion  ( - ),  the  end  of 

the  surging  motion  ( - ),  and  the  end  of  both  motions  when  the  rates  are  the  same  ( - ).  The  opaque  regions  of  the 

plots  are  regions  where  the  model  is  not  designed  to  predict  the  unsteady  force  exerted  on  the  plate. 


Along  the  diagonal  of  Figure  12  where  the  pitch-rate  matches  the  surge-rate  (cp  =  Cg),  the  model  captures 
the  increase  of  the  force  but  slightly  under  predicts  the  peak  normal  force  at  the  end  of  the  motion.  For 
the  cases  where  the  pitch-rate  is  greater  than  the  surge-rate  (cp  <  Cg),  the  model  predicts  the  similar  trend 
of  how  the  force  increase  but  under  predicts  the  peak  force  at  the  end  of  the  motion.  Finally,  for  the  cases 
in  which  the  surge-rate  is  greater  than  the  pitch-rate  (cp  >  Cg),  the  model  shows  more  deviation  from  the 
actual  forces  determined  by  the  DNS.  When  the  plates  translational  acceleration  is  faster  than  the  rotational 
acceleration,  a  quick  drop  in  the  normal  force  is  observed  before  the  motion  has  completed.  This  sudden 
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drop  is  due  to  the  development  of  the  LEV  over  the  plate.  These  trends  are  somewhat  predicted  as  the 
present  force  model  is  originally  derived  for  cases  with  the  pitch  and  surge  rates  matching.  Even  with  the 
observed  difference  off  the  diagonal  of  Figure  12,  we  note  that  the  current  predictions  are  reasonable  as  the 
motion  amplitudes  are  large  with  high  frequency  components  in  most  cases. 


V.  Comparison  with  Wagner’s  Theory 

In  order  to  assess  the  validity  for  the  force  model  developed  in  the  present  study,  the  forces  predicted 
in  the  previous  sections  are  now  compared  to  Wagner’s  inviscid  theory  for  unsteady  airfoils  at  small  angles 
of  attack.  Chen  et  al.^^  and  Pullin  and  Wang^*’  have  used  this  model  as  a  validation  for  accelerating  flat 
plates  and  ellipses’,  respectively.  For  surging  motions  with  constant  acceleration,  Wagner’s  theory  predicts 
the  forces  in  the  x  and  y  directions  for  small  s  and  a  as: 

F^,  =  |Knt/(t))  (14) 

Fy,  =  (t)sin(a)  +  (15) 


In  the  previous  two  equations,  the  primed  coordinates,  {x' ,y'),  are  in  the  direction  of  the  maneuvering 
wing,  rotated  from  the  body  fitted  coordinates,  {x,y),  and  mC  is  the  added  mass  tensor  in  these  rotated 
coordinates  (m'  =  RmR“^),  where  m  is  the  added  mass  tensor  in  the  body  fitted  coordinates  and  R  is  the 
rotation  matrix  between  the  body  fitted  coordinates  and  the  primed  coordinates,  both  given  by 


and 


(cos(a)  — sin(Q:) 
sin(Q:)  cos(a) 


(16) 


(17) 


For  the  pure  surging  motions,  the  terms  dependent  on  a  in  the  added  mass  tensor  are  independent  of 
time  (d  =  0),  therefore  Wagner’s  theory  predicts  the  forces  in  the  {x',y')  direction  as 


^  TT  2  ■  2,  \dU 

F^'  =  2^" 


(18) 


Fy,= 


TT  ttQ.  •  /  \  2  f  \  •  /  \ 

— sin(a)  +2^*^  cos(a)sm(a)  — 


(19) 


For  the  combined  maneuvers  that  include  both  surging  and  pitching  motions,  the  angle  of  attack  is  time 
dependent  (a  =  a(t)),  thus  the  derivatives  of  the  added  mass  tensor  must  be  evaluated,  and  the  resulting 
forces  in  the  {x'  ,y')  direction  are  given  by 


Px'  =  -7:  pc 


da 


dU 


2cos(a)sin(a)  — t/(t)  +  sin^(a)  — 


dt 


dt 


Fy'  =  —pcU^sin(a)  +  —pc^ 


da 


dU 


(cos"^(q;)  —  sim(a))  — C/(t)  +  cos(Q:)sin(Q:)  — 


dt 


dt 


(20) 

(21) 


In  the  analysis  to  follow,  the  forces  predicted  by  Wagner’s  inviscid  theory  are  presented  as  the  normal  force 
to  be  consistent  with  the  previous  analysis. 

We  first  consider  the  comparison  between  the  numerical  simulations,  the  present  force  modeling  results, 
and  the  theory  of  Wagner  for  the  pure  surging  case,  as  seen  in  Figure  13.  For  this  comparison,  we  consider 
the  cases  of  constant  acceleration  over  1,2,4,  and  6  chord  lengths  of  travel  Cs  =  1,2,4,  6  at  an  angle  of 
a  =  45°.  For  the  two  fastest  cases  considered  (cs  =  1  and  2),  the  present  model  predicts  the  normal  force 
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very  well  throughout  the  acceleration  as  well  as  the  peak  normal  force  at  the  end  of  the  motion.  For  the  two 
slower  cases,  the  model  begins  to  deviate  towards  the  end  of  the  acceleration,  which  could  be  due  to  the  slow 
dynamics  associated  with  the  development  of  the  LEV  over  the  wing.  It  can  be  seen  that  Wagner’s  model 
under  predicts  the  normal  force  for  all  of  the  motions  considered,  which  was  expected  due  to  the  absence  of 
the  LEV  formation  in  the  model. 
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Figure  13:  Comparison  of  the  computational  results,  the  present  force  model,  and  the  theory  of  Wagner  for  the  pure  surging 
motion  of  constant  acceleration  over  1,  2, 4,  and  6  chord  lengths  of  travel  at  Re  =  500  and  a  =  45°. 


Next,  we  assess  the  results  for  the  pure  surging  case  with  a  variation  in  angle  of  attack  by  comparing  the 
simulations,  present  force  model,  and  the  theory  of  Wagner,  as  seen  in  Figure  14.  The  constant  acceleration 
for  the  present  comparison  occurs  over  1  chord  length  of  travel  (c^  =  1)  and  the  variation  in  angle  of  attack 
is  between  a  =  15°, 30°, 45°  and  60°.  The  present  force  model  predicts  the  normal  force  throughout  the 
accelerated  portion  for  this  fast  motion  very  well  over  all  angles  of  attack  considered.  After  the  acceleration 
has  ended,  there  is  no  drop  in  the  force  that  is  predicted  with  the  present  model  when  compared  to  Wagner’s 
model.  It  can  be  seen  that  for  the  smallest  angle  considered  (a  =  15°),  Wagner’s  theory  predicts  the  normal 
force  very  well  though  the  constant  accelerated  motion  due  to  the  satisfying  the  small  angle  assumption  over 
short  travel  distances,  in  which  a  large  LEV  is  not  observed. 
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Figure  14:  Comparison  of  the  computational  results,  the  present  force  model,  and  the  theory  of  Wagner  for  the  pure  surging 
motion  of  constant  acceleration  over  1  chord  length  of  travel  and  a  variation  of  angle  of  attack  for  a  =  15°,  30°,  45°  and  60°. 


We  now  compare  the  present  force  model  for  the  combined  pitch  and  surge  motion  with  the  simulations 
and  the  theory  of  Wagner,  as  seen  in  Figure  15.  The  combined  motion  for  the  present  comparison  occurs 
over  1,2,4,  and  6  chord  lengths  of  travel  (cg  =  Cp  =  1,2,4,  and  6,  respectively)  and  the  flat  plate  pitches 
from  a  =  0°  to  ttmax  =  45°.  For  these  combined  motions,  the  present  model  predicts  the  development  of 
the  normal  force  on  the  flat  plate  for  these  complex  maneuvers  relatively  well  for  the  two  fastest  motions, 
but  under  predicts  the  peak  force  at  the  end  of  the  motion.  It  can  also  be  seen  that  after  the  motion  has 
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finished,  the  normal  force  drops  has  been  previously  observed  for  Wagner’s  model  for  the  pure  surging  cases. 
Although  the  present  model  does  under  predict  the  normal  force  for  these  motions,  it  still  reasonably  predicts 
the  forces  compared  to  the  Wagner’s  model. 
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Figure  15:  Comparison  of  the  computational  results,  the  present  force  model,  and  the  theory  of  Wagner  for  the  combined 
pitch-surge  motion  occurring  over  1,  2,  4,  and  6  chord  lengths  of  travel  (cs  =  Cp  =  1,  2, 4,  and  6,  respectively)  at  Re  =  500.  The 
pitching  motion  for  the  present  results  is  from  a  =  0°  to  Omax  =  45°. 


VI.  Discussion  on  the  Force  Model 


Let  us  further  comment  on  the  final  form  of  the  present  unsteady  nonlinear  force  model.  We  note  that 
the  model  is  comprised  of  mainly  three  components,  namely  the  parts  that  capture  surging,  pitching,  and 
incident  angle  effects,  as  shown  below: 

Fat  «  ^^sina  +  C2S^’^s  1  +  dids^  +  d2cis^  .  (22) 

V2  ' - . - "  ' - . - " 

-v  surging  pitching 

incidence  angle 

We  can  non-dimensionalize  the  above  normal  force  with  instantaneous  velocity  s 

Cn  =  ,  «  •\/2sina[cis^^“^  +  C2S^’^“^s]  +  2did  +  2^20,  (23) 

where  we  have  taken  poo  and  c  to  be  unity  without  the  loss  of  generality. 

Next,  let  us  compare  our  model  with  other  linear  aerodynamic  force  models.  For  a  wing  at  fixed  low 
angle  of  attack  with  steady  motion  (in  absence  of  pitching),  we  find  that  sin  a  «  a  and  the  above  expression 
becomes 

Cn  «  y2cis^i-2a  =  4.21s°  “a  «  4.21a  =  0.67(27r)a,  (24) 

which  tells  us  that  the  present  model  estimates  the  steady  lift  slope  Og  =  dCi/da  to  be  0.67%  of  the  potential 
value  (27r).  This  is  in  agreement  with  the  numerical  finding  from  Taira  and  Colonius"^  that  reports  the  lift 
close  at  low  a  to  be  og  =  0.73(27r).  The  lower  lift  slope  compared  to  the  potential  value  is  due  to  the  viscous 
effects  at  the  low  Reynolds  number  selected  here.  For  higher  Reynolds  number,  Ci  can  be  replaced  with  a 
larger  value  to  take  Reynolds  number  effect  into  consideration  (i.e.  ci  =  ci(i?e)). 

We  also  note  that  for  small  s  and  a,  the  pitch-surge  model  reduces  to 

Cn  ^  V2ci{ks)a  +  2di{kp)a  +  2d2a,  (25) 


where  dependence  on  pitch  and  surge-based  reduced  frequencies  kp  and  kg,  respectively,  are  shown  if  present. 
It  can  be  noticed  that  this  form  of  the  equation  shares  similarity  with  Theodorsen’s  modeb 


Cl  •  Theodorsen  —  27rC{k)a  + 


TT 

4 


37r  c 


C{k) 


TT  C 


(26) 
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where  C{k)  is  the  Theodorsen’s  function.^ 

In  upcoming  studies,  we  will  incorporate  the  effect  of  aspect  ratio  using  the  a  relationship  similar  to  the 
one  proposed  by  Helmbold^^ 


^  Cl, 2D 

C^L,3D  —  --  -- 

^Jl  +  {2/ARf  +  2/AR 


(27) 


into  the  unsteady  model  developed  for  the  two-dimensional  wing. 


VI I .  S  ummar  y 

An  nonlinear  aerodynamics  force  model  is  developed  for  wings  undergoing  large-amplitude  unsteady 
maneuvers.  The  motion  considered  includes  surging,  pitching,  and  the  combination  of  the  two  for  low 
Reynolds  number  two-dimensional  flow.  The  proposed  model  is  able  to  accurately  predict  normal  force 
on  the  wing  based  on  two-dimensional  DNS  results  over  a  range  of  motion  profiles  and  angles  of  attack. 
One  advantage  of  the  present  model  is  that  the  input  is  only  the  wing  kinematics  and  does  not  require 
any  knowledge  of  vortex  dynamics.  For  this  reason,  the  normal  force  can  be  predicted  with  practically 
no  computational  expense.  What  remains  to  be  completed  with  this  research  effort  is  to  further  test  the 
proposed  model  for  a  wider  variety  of  wing  motion  and  higher  Reynolds  number  flows.  Since  knowledge 
from  past  research^^’  suggests  that  force  on  the  wing  is  fairly  Reynolds  number  independent  for  flow  below 
transition  {Re  <  10^),  we  expect  the  present  force  model  to  be  applicable  to  cases  with  higher  Reynolds 
numbers. 
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